Surface emitter
In [1]:
Copied!
# Useful for debugging
%load_ext autoreload
%autoreload 2
%config InlineBackend.figure_format = 'retina'
import os
import yaml
import numpy as np
# Useful for debugging
%load_ext autoreload
%autoreload 2
%config InlineBackend.figure_format = 'retina'
import os
import yaml
import numpy as np
In [2]:
Copied!
from distgen import Generator
from distgen import Generator
In [3]:
Copied!
input_yaml = """
n_particle: 100000
species: electron
start:
type: free
random:
type: hammersley
total_charge:
units: pC
value: 10
r_dist:
max_r:
units: mm
value: 1
type: radial_uniform
p_dist:
sigma_p:
value: 1
units: eV/c
avg_p:
value: 100
units: eV/c
type: gaussian
"""
input_yaml = """
n_particle: 100000
species: electron
start:
type: free
random:
type: hammersley
total_charge:
units: pC
value: 10
r_dist:
max_r:
units: mm
value: 1
type: radial_uniform
p_dist:
sigma_p:
value: 1
units: eV/c
avg_p:
value: 100
units: eV/c
type: gaussian
"""
In [4]:
Copied!
gen = Generator(input_yaml, verbose=True)
gen = Generator(input_yaml, verbose=True)
In [5]:
Copied!
P = gen.run()
P = gen.run()
Distribution format: None
Warning: no output file specified, defaulting to "None".
Output file: None
Creating beam distribution....
Beam starting from: free
Total charge: 10 pC.
Number of macroparticles: 100000.
Assuming cylindrical symmetry...
r distribution: radial uniform
min_r = 0 mm, max_r = 1 mm
theta distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
p distribution: Gaussian
avg_p = 100 eV/c, sigma_p = 1.000 eV/c
azimuthal angle distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
polar angle distribution: uniform phi
min_phi = 0 rad, max_phi = 3.14159 rad
Shifting avg_x = -4.85224E-08 mm -> 0 mm
Scaling sigma_x = 0.499991 mm -> 0.5 mm
Shifting avg_y = 1.44675E-06 mm -> 0 mm
Scaling sigma_y = 0.49999 mm -> 0.5 mm
Shifting avg_px = -3.3074E-05 eV/c -> -6.12323E-15 eV/c
Scaling sigma_px = 57.7382 eV/c -> 57.7379 eV/c
Shifting avg_py = 1.54491E-06 eV/c -> 0 eV/c
Scaling sigma_py = 57.7382 eV/c -> 57.7379 eV/c
Shifting avg_pz = -0.000118434 eV/c -> 6.12323E-15 eV/c
Scaling sigma_pz = 57.7372 eV/c -> 57.7379 eV/c
...done. Time Elapsed: 104.218 ms.
Created particles in .particles:
ParticleGroup with 100000 particles with total charge 1.0000000000000003e-11 C
In [6]:
Copied!
P.plot('p')
P.plot('p')
In [7]:
Copied!
P.plot('x', 'px')
P.plot('x', 'px')
In [8]:
Copied!
P.plot('pz')
P.plot('pz')
In [9]:
Copied!
P.plot('px')
P.plot('px')
In [10]:
Copied!
P.plot('py')
P.plot('py')
In [11]:
Copied!
from matplotlib import pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(projection='3d')
px, py, pz = P['px'], P['py'], P['pz']
p = np.sqrt(px**2 + py**2 + pz**2)
theta = np.arctan2(py, px)
phi = np.arccos(pz/p)
hist, pedges = np.histogram(phi, bins=100, density=True)
pcs = (pedges[1:] + pedges[:-1]) / 2
x = px/p
y = py/p
z = pz/p
ax.scatter(x[::200], y[::200], z[::200], '.');
ax.set_xlabel(r'$\hat{p}_x$')
ax.set_ylabel(r'$\hat{p}_y$')
ax.set_zlabel(r'$\hat{p}_z$')
from matplotlib import pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(projection='3d')
px, py, pz = P['px'], P['py'], P['pz']
p = np.sqrt(px**2 + py**2 + pz**2)
theta = np.arctan2(py, px)
phi = np.arccos(pz/p)
hist, pedges = np.histogram(phi, bins=100, density=True)
pcs = (pedges[1:] + pedges[:-1]) / 2
x = px/p
y = py/p
z = pz/p
ax.scatter(x[::200], y[::200], z[::200], '.');
ax.set_xlabel(r'$\hat{p}_x$')
ax.set_ylabel(r'$\hat{p}_y$')
ax.set_zlabel(r'$\hat{p}_z$')
Out[11]:
Text(0.5, 0, '$\\hat{p}_z$')
In [12]:
Copied!
theta = np.arctan2(py, px)
phi = np.arccos(pz/p)
hist, pedges = np.histogram(phi, bins=100, density=True)
pcs = (pedges[1:] + pedges[:-1]) / 2
plt.plot(pcs, hist);
theta = np.arctan2(py, px)
phi = np.arccos(pz/p)
hist, pedges = np.histogram(phi, bins=100, density=True)
pcs = (pedges[1:] + pedges[:-1]) / 2
plt.plot(pcs, hist);
In [13]:
Copied!
rns = np.random.random(100000)
rns = np.random.random(100000)
In [14]:
Copied!
pa = 0
pb = np.pi
Ca = np.cos(pa)
Cb = np.cos(pb)
pa = 0
pb = np.pi
Ca = np.cos(pa)
Cb = np.cos(pb)
In [15]:
Copied!
rho = 1/(Ca-Cb)
rho = 1/(Ca-Cb)
In [16]:
Copied!
phis = np.linspace(pa, pb, 1000)
phis = np.linspace(pa, pb, 1000)
In [17]:
Copied!
plt.plot(phis, rho*np.sin(phis));
plt.plot(phis, rho*np.sin(phis));
In [18]:
Copied!
np.trapz(rho*np.sin(phis), phis)
np.trapz(rho*np.sin(phis), phis)
/tmp/ipykernel_2251/4275813628.py:1: DeprecationWarning: `trapz` is deprecated. Use `trapezoid` instead, or one of the numerical integration functions in `scipy.integrate`. np.trapz(rho*np.sin(phis), phis)
Out[18]:
np.float64(0.999999175885426)
In [19]:
Copied!
cdf = (Ca - np.cos(phis))*rho
cdf = (Ca - np.cos(phis))*rho
In [20]:
Copied!
from scipy.integrate import cumulative_trapezoid as cumtrapz
from scipy.integrate import cumulative_trapezoid as cumtrapz
In [21]:
Copied!
plt.plot(phis, cdf, phis, cumtrapz(rho*np.sin(phis), phis, initial=0));
plt.plot(phis, cdf, phis, cumtrapz(rho*np.sin(phis), phis, initial=0));
In [22]:
Copied!
ps = np.arccos( Ca - rns*(Ca-Cb) )
ps = np.arccos( Ca - rns*(Ca-Cb) )
In [23]:
Copied!
hist, pedges = np.histogram(ps, bins=30, density=True)
hist, pedges = np.histogram(ps, bins=30, density=True)
In [24]:
Copied!
pcs = (pedges[1:] + pedges[:-1]) / 2
pcs = (pedges[1:] + pedges[:-1]) / 2
In [25]:
Copied!
plt.plot(pcs, hist, phis, rho*np.sin(phis));
plt.plot(pcs, hist, phis, rho*np.sin(phis));
In [26]:
Copied!
gen = Generator('data/beer.can.in.yaml', verbose=1)
#print(gen)
gen = Generator('data/beer.can.in.yaml', verbose=1)
#print(gen)
In [27]:
Copied!
P1 = gen.run()
P1 = gen.run()
Distribution format: gpt
Output file: beer.can.out.txt
Creating beam distribution....
Beam starting from: cathode
Total charge: 10 pC.
Number of macroparticles: 200000.
Assuming cylindrical symmetry...
r distribution: radial uniform
min_r = 0 mm, max_r = 2 mm
theta distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
t distribution: uniform
min_t = -2 ps, max_t = 2 ps, avg_t = 0 ps, sigma_t: 1.1547 ps
px distribution: Gaussian
avg_px = 0 eV/c, sigma_px = 276.857 eV/c
py distribution: Gaussian
avg_py = 0 eV/c, sigma_py = 276.857 eV/c
pz distribution: Gaussian
avg_pz = 0 eV/c, sigma_pz = 276.857 eV/c
Shifting avg_x = -2.10539E-05 mm -> 0 mm
Scaling sigma_x = 0.999982 mm -> 1 mm
Shifting avg_y = -1.97648E-06 mm -> 0 mm
Scaling sigma_y = 0.999998 mm -> 1 mm
Shifting avg_px = -0.0212516 eV/c -> 0 eV/c
Scaling sigma_px = 276.849 eV/c -> 276.857 eV/c
Shifting avg_py = -0.0259627 eV/c -> 0 eV/c
Scaling sigma_py = 276.844 eV/c -> 276.857 eV/c
Shifting avg_pz = -0.0372953 eV/c -> 0 eV/c
Scaling sigma_pz = 276.843 eV/c -> 276.857 eV/c
Shifting avg_t = -4.90022E-05 ps -> 0 ps
Scaling sigma_t = 1.1547 ps -> 1.1547 ps
Cathode start: fixing pz momenta to forward hemisphere
avg_pz -> 220.904 eV/c, sigma_pz -> 166.887 eV/c
...done. Time Elapsed: 191.009 ms.
Created particles in .particles: ParticleGroup with 200000 particles with total charge 1.0000000000000003e-11 C
In [28]:
Copied!
P1.plot('x', 'px')
P1.plot('x', 'px')
In [29]:
Copied!
P1.plot('p')
P1.plot('p')
In [30]:
Copied!
P1.plot('pz')
P1.plot('pz')
In [31]:
Copied!
P1['sigma_px'], P1['sigma_py'], P1['sigma_pz']
P1['sigma_px'], P1['sigma_py'], P1['sigma_pz']
Out[31]:
(np.float64(276.8570795554991), np.float64(276.8570795554991), np.float64(166.88683722575215))
In [32]:
Copied!
gen = Generator('data/maxwell_boltzmann.beer.can.in.yaml', verbose=1)
gen = Generator('data/maxwell_boltzmann.beer.can.in.yaml', verbose=1)
In [33]:
Copied!
P2 = gen.run()
P2 = gen.run()
Distribution format: gpt
Output file: beer.can.out.txt
Creating beam distribution....
Beam starting from: cathode
Total charge: 10 pC.
Number of macroparticles: 200000.
Assuming cylindrical symmetry...
r distribution: radial uniform
min_r = 0 mm, max_r = 2 mm
theta distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
p distribution: Maxwell-Boltzmann
p scale = 276.857 eV/c
azimuthal angle distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
polar angle distribution: uniform phi
min_phi = 0 rad, max_phi = 3.14159 rad
t distribution: uniform
min_t = -2 ps, max_t = 2 ps, avg_t = 0 ps, sigma_t: 1.1547 ps
Shifting avg_x = 1.63849E-05 mm -> 0 mm
Scaling sigma_x = 0.999984 mm -> 1 mm
Shifting avg_y = -2.77139E-06 mm -> 0 mm
Scaling sigma_y = 0.999996 mm -> 1 mm
Shifting avg_px = -0.000954618 eV/c -> -2.70524E-14 eV/c
Scaling sigma_px = 276.853 eV/c -> 276.857 eV/rad⁰⋅⁵/c
Shifting avg_py = 0.000856608 eV/c -> 0 eV/c
Scaling sigma_py = 276.848 eV/c -> 276.857 eV/rad⁰⋅⁵/c
Shifting avg_pz = -0.0142766 eV/c -> 2.70524E-14 eV/c
Scaling sigma_pz = 276.833 eV/c -> 276.857 eV/rad⁰⋅⁵/c
Shifting avg_t = -4.90022E-05 ps -> 0 ps
Scaling sigma_t = 1.1547 ps -> 1.1547 ps
Cathode start: fixing pz momenta to forward hemisphere
avg_pz -> 220.909 eV/c, sigma_pz -> 166.881 eV/c
...done. Time Elapsed: 187.787 ms.
Created particles in .particles: ParticleGroup with 200000 particles with total charge 1.0000000000000003e-11 C
In [34]:
Copied!
P2['sigma_px'], P2['sigma_py'], P2['sigma_pz']
P2['sigma_px'], P2['sigma_py'], P2['sigma_pz']
Out[34]:
(np.float64(276.85707955549907), np.float64(276.85707955549907), np.float64(166.88083659251544))
In [35]:
Copied!
gen = Generator('data/maxwell_boltzmann_KE.beer.can.in.yaml', verbose=1)
gen = Generator('data/maxwell_boltzmann_KE.beer.can.in.yaml', verbose=1)
In [36]:
Copied!
P3 = gen.run()
P3 = gen.run()
Distribution format: gpt
Output file: beer.can.out.txt
Creating beam distribution....
Beam starting from: cathode
Total charge: 10 pC.
Number of macroparticles: 200000.
Assuming cylindrical symmetry...
r distribution: radial uniform
min_r = 0 mm, max_r = 2 mm
theta distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
KE distribution: Maxwell-Boltzmann Energy
kT = 150 meV
azimuthal angle distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
polar angle distribution: uniform phi
min_phi = 0 rad, max_phi = 3.14159 rad
t distribution: uniform
min_t = -2 ps, max_t = 2 ps, avg_t = 0 ps, sigma_t: 1.1547 ps
Shifting avg_x = 1.63849E-05 mm -> 0 mm
Scaling sigma_x = 0.999984 mm -> 1 mm
Shifting avg_y = -2.77139E-06 mm -> 0 mm
Scaling sigma_y = 0.999996 mm -> 1 mm
Shifting avg_px = -1.75487E-09 meV·s/m -> -9.02058E-20 meV·s/m
Scaling sigma_px = 0.000922998 meV·s/m -> 0.000922972 meV·s/m
Shifting avg_py = 1.21968E-09 meV·s/m -> 0 meV·s/m
Scaling sigma_py = 0.000922978 meV·s/m -> 0.000922972 meV·s/m
Shifting avg_pz = -4.48933E-08 meV·s/m -> 9.02058E-20 meV·s/m
Scaling sigma_pz = 0.00092294 meV·s/m -> 0.000922972 meV·s/m
Shifting avg_t = -4.90022E-05 ps -> 0 ps
Scaling sigma_t = 1.1547 ps -> 1.1547 ps
Cathode start: fixing pz momenta to forward hemisphere
avg_pz -> 0.0007366 meV·s/m, sigma_pz -> 0.000556145 meV·s/m
...done. Time Elapsed: 190.986 ms.
Created particles in .particles: ParticleGroup with 200000 particles with total charge 1.0000000000000003e-11 C
In [37]:
Copied!
P3.plot('kinetic_energy')
P3.plot('kinetic_energy')
In [38]:
Copied!
P3.plot('x', 'px')
P3.plot('x', 'px')
In [39]:
Copied!
P3['sigma_px'], P3['sigma_py'], P3['sigma_pz']
P3['sigma_px'], P3['sigma_py'], P3['sigma_pz']
Out[39]:
(np.float64(276.70000872169317), np.float64(276.70000872169317), np.float64(166.7279863589753))
In [40]:
Copied!
fig = plt.figure()
ax = fig.add_subplot(projection='3d')
px, py, pz = P3['px'], P3['py'], P3['pz']
p = np.sqrt(px**2 + py**2 + pz**2)
theta = np.arctan2(py, px)
phi = np.arccos(pz/p)
hist, pedges = np.histogram(phi, bins=100, density=True)
pcs = (pedges[1:] + pedges[:-1]) / 2
x = px/p
y = py/p
z = pz/p
ax.scatter(x[::500], y[::500], z[::500], '.');
ax.set_xlabel(r'$\hat{p}_x$')
ax.set_ylabel(r'$\hat{p}_y$')
ax.set_zlabel(r'$\hat{p}_z$')
fig = plt.figure()
ax = fig.add_subplot(projection='3d')
px, py, pz = P3['px'], P3['py'], P3['pz']
p = np.sqrt(px**2 + py**2 + pz**2)
theta = np.arctan2(py, px)
phi = np.arccos(pz/p)
hist, pedges = np.histogram(phi, bins=100, density=True)
pcs = (pedges[1:] + pedges[:-1]) / 2
x = px/p
y = py/p
z = pz/p
ax.scatter(x[::500], y[::500], z[::500], '.');
ax.set_xlabel(r'$\hat{p}_x$')
ax.set_ylabel(r'$\hat{p}_y$')
ax.set_zlabel(r'$\hat{p}_z$')
Out[40]:
Text(0.5, 0, '$\\hat{p}_z$')
In [41]:
Copied!
KE = np.linspace(0, 150, 10000)
E0=0
A1 = 0.8,
m1 = 8
A2 = 0.1
m2 = 90
PKE = A1*np.exp(-np.abs(KE-E0)/m1) + A2*np.exp(-np.abs(KE-E0)/m2)
PKE = PKE/np.trapz(PKE, KE) # Numerically intergate to normalize
plt.plot(KE, PKE);
plt.xlabel('KE (eV)');
plt.ylabel('$\\rho(KE)$ (1/eV)');
KE = np.linspace(0, 150, 10000)
E0=0
A1 = 0.8,
m1 = 8
A2 = 0.1
m2 = 90
PKE = A1*np.exp(-np.abs(KE-E0)/m1) + A2*np.exp(-np.abs(KE-E0)/m2)
PKE = PKE/np.trapz(PKE, KE) # Numerically intergate to normalize
plt.plot(KE, PKE);
plt.xlabel('KE (eV)');
plt.ylabel('$\\rho(KE)$ (1/eV)');
/tmp/ipykernel_2251/1250638860.py:12: DeprecationWarning: `trapz` is deprecated. Use `trapezoid` instead, or one of the numerical integration functions in `scipy.integrate`. PKE = PKE/np.trapz(PKE, KE) # Numerically intergate to normalize
In [42]:
Copied!
dat = np.zeros((len(KE),2))
dat[:,0], dat[:,1] = KE, PKE
np.savetxt('KEdist.txt', dat, header="KE PKE", comments='')
dat = np.zeros((len(KE),2))
dat[:,0], dat[:,1] = KE, PKE
np.savetxt('KEdist.txt', dat, header="KE PKE", comments='')
In [43]:
Copied!
q = 1000*1.60217663e-19
print(q)
input_yaml = """
n_particle: 1000
species: electron
start:
type: cathode
random:
type: hammersley
total_charge:
units: C
value: 1.60217663e-17
r_dist:
max_r:
units: mm
value: 14.6
type: radial_uniform
KE_dist:
file: KEdist.txt
units: eV
type: file1d
transforms:
sx:
avg_x:
units: millimeter
value: -50
type: set_avg x
sy:
avg_y:
units: millimeter
value: 50
type: set_avg y
sz:
avg_z:
units: millimeter
value: 65
type: set_avg z
output:
file: None
"""
q = 1000*1.60217663e-19
print(q)
input_yaml = """
n_particle: 1000
species: electron
start:
type: cathode
random:
type: hammersley
total_charge:
units: C
value: 1.60217663e-17
r_dist:
max_r:
units: mm
value: 14.6
type: radial_uniform
KE_dist:
file: KEdist.txt
units: eV
type: file1d
transforms:
sx:
avg_x:
units: millimeter
value: -50
type: set_avg x
sy:
avg_y:
units: millimeter
value: 50
type: set_avg y
sz:
avg_z:
units: millimeter
value: 65
type: set_avg z
output:
file: None
"""
1.6021766299999998e-16
In [44]:
Copied!
import yaml
inputs = yaml.safe_load(input_yaml)
n = 1000
inputs['n_particle'] = n
inputs['output']['file'] = 'test_ion_writer.ion'
inputs['output']['type'] = 'simion'
inputs['total_charge']['value'] = n*1.60217663e-19
inputs['transforms']['sx']['avg_x']['value'] = -50
inputs['transforms']['sy']['avg_y']['value'] = 50
inputs['transforms']['sz']['avg_z']['value'] = 55
import yaml
inputs = yaml.safe_load(input_yaml)
n = 1000
inputs['n_particle'] = n
inputs['output']['file'] = 'test_ion_writer.ion'
inputs['output']['type'] = 'simion'
inputs['total_charge']['value'] = n*1.60217663e-19
inputs['transforms']['sx']['avg_x']['value'] = -50
inputs['transforms']['sy']['avg_y']['value'] = 50
inputs['transforms']['sz']['avg_z']['value'] = 55
In [45]:
Copied!
gen = Generator(inputs, verbose=1)
gen = Generator(inputs, verbose=1)
In [46]:
Copied!
P = gen.run()
P = gen.run()
Distribution format: simion
Output file: /home/runner/work/distgen/distgen/docs/examples/test_ion_writer.ion
Creating beam distribution....
Beam starting from: cathode
Total charge: 1.60218E-16 C.
Number of macroparticles: 1000.
Assuming cylindrical symmetry...
r distribution: radial uniform
min_r = 0 mm, max_r = 14.6 mm
theta distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
KE distribution: KE-distribution file: "/home/runner/work/distgen/distgen/docs/examples/KEdist.txt"
azimuthal angle distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
polar angle distribution: uniform phi
min_phi = 0 rad, max_phi = 3.14159 rad
Shifting avg_x = 0.00564031 mm -> 0 mm
Scaling sigma_x = 7.29551 mm -> 7.3 mm
Shifting avg_y = 0.00120927 mm -> 0 mm
Scaling sigma_y = 7.29741 mm -> 7.3 mm
Shifting avg_px = -6.09774E-09 eV·s/m -> -9.89015E-22 eV·s/m
Scaling sigma_px = 1.11822E-05 eV·s/m -> 1.11773E-05 eV·s/m
Shifting avg_py = -2.50205E-09 eV·s/m -> 0 eV·s/m
Scaling sigma_py = 1.11987E-05 eV·s/m -> 1.11773E-05 eV·s/m
Shifting avg_pz = -2.52855E-08 eV·s/m -> 9.89015E-22 eV·s/m
Scaling sigma_pz = 1.11509E-05 eV·s/m -> 1.11773E-05 eV·s/m
Cathode start: fixing pz momenta to forward hemisphere
avg_pz -> 8.07962E-06 eV·s/m, sigma_pz -> 7.72348E-06 eV·s/m
Applying user supplied transform: "sx" = set_avg x...
Setting avg_x -> -50 mm.
Applying user supplied transform: "sy" = set_avg y...
Setting avg_y -> 50 mm.
Applying user supplied transform: "sz" = set_avg z...
Setting avg_z -> 55 mm.
...done. Time Elapsed: 26.2625 ms.
Created particles in .particles:
ParticleGroup with 1000 particles with total charge 1.6021766300000003e-16 C
In [47]:
Copied!
P.plot('x', 'y', figsize=(5,5))
P.plot('x', 'y', figsize=(5,5))
In [48]:
Copied!
P.plot('x', 'px', figsize=(5,5))
P.plot('x', 'px', figsize=(5,5))
In [49]:
Copied!
P.plot('kinetic_energy')
P.plot('kinetic_energy')
In [50]:
Copied!
from distgen.writers import writer
writer('simion', gen.beam(), 'test_ion_writer.ion', params={'color':0})
from distgen.writers import writer
writer('simion', gen.beam(), 'test_ion_writer.ion', params={'color':0})
Distribution format: simion
Output file: /home/runner/work/distgen/distgen/docs/examples/test_ion_writer.ion
Creating beam distribution....
Beam starting from: cathode
Total charge: 1.60218E-16 C.
Number of macroparticles: 1000.
Assuming cylindrical symmetry...
r distribution: radial uniform
min_r = 0 mm, max_r = 14.6 mm
theta distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
KE distribution: KE-distribution file: "/home/runner/work/distgen/distgen/docs/examples/KEdist.txt"
azimuthal angle distribution: uniform theta
min_theta = 0 rad, max_theta = 6.28319 rad
polar angle distribution: uniform phi
min_phi = 0 rad, max_phi = 3.14159 rad
Shifting avg_x = 0.00564031 mm -> 0 mm
Scaling sigma_x = 7.29551 mm -> 7.3 mm
Shifting avg_y = 0.00120927 mm -> 0 mm
Scaling sigma_y = 7.29741 mm -> 7.3 mm
Shifting avg_px = -6.09774E-09 eV·s/m -> -9.89015E-22 eV·s/m
Scaling sigma_px = 1.11822E-05 eV·s/m -> 1.11773E-05 eV·s/m
Shifting avg_py = -2.50205E-09 eV·s/m -> 0 eV·s/m
Scaling sigma_py = 1.11987E-05 eV·s/m -> 1.11773E-05 eV·s/m
Shifting avg_pz = -2.52855E-08 eV·s/m -> 9.89015E-22 eV·s/m
Scaling sigma_pz = 1.11509E-05 eV·s/m -> 1.11773E-05 eV·s/m
Cathode start: fixing pz momenta to forward hemisphere
avg_pz -> 8.07962E-06 eV·s/m, sigma_pz -> 7.72348E-06 eV·s/m
Applying user supplied transform: "sx" = set_avg x...
Setting avg_x -> -50 mm.
Applying user supplied transform: "sy" = set_avg y...
Setting avg_y -> 50 mm.
Applying user supplied transform: "sz" = set_avg z...
Setting avg_z -> 55 mm.
...done. Time Elapsed: 26.4134 ms.
In [51]:
Copied!
from scipy.constants import physical_constants
mc2 = 1e6 * physical_constants['electron mass energy equivalent in MeV'][0]
e_= physical_constants['elementary charge'][0]
me = physical_constants['electron mass in u'][0]
def particle_group_to_SIMION(P, filename, color=0):
header=';0'
simion_params= ['TOB', 'MASS', 'CHARGE', 'X', 'Y', 'Z', 'AZ', 'EL', 'KE', 'CWF', 'COLOR']
simion_units = {'TOB':'usec', 'MASS':'amu', 'CHARGE':'e', 'X':'mm', 'Y':'mm', 'Z':'mm', 'AZ':'deg', 'EL':'deg', 'CWF':'', 'COLOR':''}
data = np.zeros( (len(P), len(simion_params)) )
data[:, simion_params.index('TOB')] = P.t*1e6 # [P.t] = sec, convert to usec
if(P.species == 'electron'):
data[:, simion_params.index('MASS')] = np.full(len(P), me)
data[:, simion_params.index('CHARGE')] = np.full(len(P), -1)
else:
raise ValueError(f'Species {P.species} is not supported')
data[:, simion_params.index('X')] = P.z*1e3
data[:, simion_params.index('Y')] = P.y*1e3
data[:, simion_params.index('Z')] = -P.x*1e3
px = P.pz
py = P.py
pz = -P.px
data[:, simion_params.index('KE')] = P.kinetic_energy # [eV]
data[:, simion_params.index('AZ')] = np.arctan2(-pz, px) * (180/np.pi) # [deg]
data[:, simion_params.index('EL')] = np.arctan2(py, np.sqrt(px**2 + pz**2) ) * (180/np.pi) # [deg]
data[:, simion_params.index('CWF')] = P.weight/e_ # Charge Weighting Factor, derive from particle group weights
data[:, simion_params.index('COLOR')] = np.full(len(P), color)
#fname, X, fmt='%.18e', delimiter=' '
np.savetxt(filename, data, delimiter=',', header=header, comments='', fmt=' %.9e')
from scipy.constants import physical_constants
mc2 = 1e6 * physical_constants['electron mass energy equivalent in MeV'][0]
e_= physical_constants['elementary charge'][0]
me = physical_constants['electron mass in u'][0]
def particle_group_to_SIMION(P, filename, color=0):
header=';0'
simion_params= ['TOB', 'MASS', 'CHARGE', 'X', 'Y', 'Z', 'AZ', 'EL', 'KE', 'CWF', 'COLOR']
simion_units = {'TOB':'usec', 'MASS':'amu', 'CHARGE':'e', 'X':'mm', 'Y':'mm', 'Z':'mm', 'AZ':'deg', 'EL':'deg', 'CWF':'', 'COLOR':''}
data = np.zeros( (len(P), len(simion_params)) )
data[:, simion_params.index('TOB')] = P.t*1e6 # [P.t] = sec, convert to usec
if(P.species == 'electron'):
data[:, simion_params.index('MASS')] = np.full(len(P), me)
data[:, simion_params.index('CHARGE')] = np.full(len(P), -1)
else:
raise ValueError(f'Species {P.species} is not supported')
data[:, simion_params.index('X')] = P.z*1e3
data[:, simion_params.index('Y')] = P.y*1e3
data[:, simion_params.index('Z')] = -P.x*1e3
px = P.pz
py = P.py
pz = -P.px
data[:, simion_params.index('KE')] = P.kinetic_energy # [eV]
data[:, simion_params.index('AZ')] = np.arctan2(-pz, px) * (180/np.pi) # [deg]
data[:, simion_params.index('EL')] = np.arctan2(py, np.sqrt(px**2 + pz**2) ) * (180/np.pi) # [deg]
data[:, simion_params.index('CWF')] = P.weight/e_ # Charge Weighting Factor, derive from particle group weights
data[:, simion_params.index('COLOR')] = np.full(len(P), color)
#fname, X, fmt='%.18e', delimiter=' '
np.savetxt(filename, data, delimiter=',', header=header, comments='', fmt=' %.9e')
In [52]:
Copied!
particle_group_to_SIMION(P, 'text_ion_file.ion')
particle_group_to_SIMION(P, 'text_ion_file.ion')
In [53]:
Copied!
def read_simion_ION_file(filename):
data = np.loadtxt(filename, comments=';', delimiter=',', skiprows=1)
simion_params= ['TOB', 'MASS', 'CHARGE', 'X', 'Y', 'Z', 'AZ', 'EL', 'KE', 'CWF', 'COLOR']
return {simion_params[ii]:data[:,ii] for ii, p in enumerate(simion_params)}
def read_simion_ION_file(filename):
data = np.loadtxt(filename, comments=';', delimiter=',', skiprows=1)
simion_params= ['TOB', 'MASS', 'CHARGE', 'X', 'Y', 'Z', 'AZ', 'EL', 'KE', 'CWF', 'COLOR']
return {simion_params[ii]:data[:,ii] for ii, p in enumerate(simion_params)}
In [54]:
Copied!
ions1 = read_simion_ION_file('text_ion_file.ion')
ions2 = read_simion_ION_file('test_ion_writer.ion')
ions1 = read_simion_ION_file('text_ion_file.ion')
ions2 = read_simion_ION_file('test_ion_writer.ion')
In [55]:
Copied!
fig = plt.figure()
ax = fig.add_subplot(projection='3d')
px, py, pz = P['px'], P['py'], P['pz']
p = np.sqrt(px**2 + py**2 + pz**2)
theta = np.arctan2(py, px)
phi = np.arccos(pz/p)
hist, pedges = np.histogram(phi, bins=100, density=True)
pcs = (pedges[1:] + pedges[:-1]) / 2
x = px/p
y = py/p
z = pz/p
ax.scatter(x[::], y[::], z[::], '.');
ax.set_xlabel(r'$\hat{p}_x$')
ax.set_ylabel(r'$\hat{p}_y$')
ax.set_zlabel(r'$\hat{p}_z$')
fig = plt.figure()
ax = fig.add_subplot(projection='3d')
px, py, pz = P['px'], P['py'], P['pz']
p = np.sqrt(px**2 + py**2 + pz**2)
theta = np.arctan2(py, px)
phi = np.arccos(pz/p)
hist, pedges = np.histogram(phi, bins=100, density=True)
pcs = (pedges[1:] + pedges[:-1]) / 2
x = px/p
y = py/p
z = pz/p
ax.scatter(x[::], y[::], z[::], '.');
ax.set_xlabel(r'$\hat{p}_x$')
ax.set_ylabel(r'$\hat{p}_y$')
ax.set_zlabel(r'$\hat{p}_z$')
Out[55]:
Text(0.5, 0, '$\\hat{p}_z$')
In [56]:
Copied!
for p in ions1: print(f'{p}:', max(np.abs(ions1[p] - ions2[p])))
for p in ions1: print(f'{p}:', max(np.abs(ions1[p] - ions2[p])))
TOB: 0.0 MASS: 0.0 CHARGE: 0.0 X: 0.0 Y: 0.0 Z: 0.0 AZ: 0.0 EL: 0.0 KE: 0.0 CWF: 0.0 COLOR: 0.0
In [57]:
Copied!
plt.plot(ions1['KE'], ions1['KE']/ions2['KE'], '.');
plt.xlabel('KE (eV)');
plt.plot(ions1['KE'], ions1['KE']/ions2['KE'], '.');
plt.xlabel('KE (eV)');
In [58]:
Copied!
os.remove('text_ion_file.ion')
os.remove('KEdist.txt')
os.remove('test_ion_writer.ion')
os.remove('text_ion_file.ion')
os.remove('KEdist.txt')
os.remove('test_ion_writer.ion')
In [ ]:
Copied!
In [ ]:
Copied!